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are  presented.  The  dissipation  factor  of  a  turbulent  air  masSscan  be  measured  by  an  aircraft  and  a  weather 
radar.  Comparisons  are  made  between  precipitation  reflectivity  and  spectral  width  measurements  as  indict. "ors 
at  turbulence.  The  instrumentation  and  data  processing  procedures  are  described.  Examples  of  turbulence 
observations  made  with  a  storm-penetrating  aircraft  and  the  we  ether  radar  are  given.  The  relationship  between 
the  radar  observations  and  the  physical  properties  of  die  urbuleqffc  atmosphere  are  derived.  The  relationship 
of  radar  spectral  width  (variance)  to  turbulence  intensity  is  discussed. 
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COORDINATED  RADAR  AND  AIRCRAFT  OBSERVATIONS  OF  TURBULENCE 
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1.0  INTRODUCTION 

The  ability  to  sense  and  display  regions  of  hazardous  turbulence  aloft  is 
an  essential  prerequisite  to  the  selection  of  safe,  minimum-distance-traveled 
flight  paths.  Lincoln  Laboratory,  under  a  program  sponsored  by  the  Federal 
Aviation  Administration,  has  undertaken  to  measure  and  correlate  radar-  and 
aircraft-sensed  turbulence  in  rainstorms.*  This  paper  presents  interim  re¬ 
sults  of  these  measurements.  Results  of  the  measurements  are  expected  to  pro¬ 
vide  guidance  in  the  design  of  equipment  that  will  display  areas  of  hazardous 
turbulence  with  confidence.  Plots  of  turbulence  levels  measured  simultaneous¬ 
ly  by  a  ground  radar  and  an  aircraft  are  shown  to  illustrate  the  correlation. 
Comparisons  of  the  turbulence  levels  with  the  rain  reflectivity  are  also  made. 
2.0  BACKGROUND 

Precipitation  reflectivity  and  spectral  width  as  measured  by  a  ground- 
based  weather  radar  have  been  used  as  turbulence  indicators.  The  current  re¬ 
search  focuses  on  an  investigation  of  the  use  of  spectral  width  as  a  reliable 
Indicator  of  turbulence  in  precipitation. 

The  dissipation  factor,  e,  of  the  turbulent  air  mass  can  be  measured  by 

JTS*  - 

both  the  aircraft  and  the  radar  and  is  independent  of  type  and  speed  of  the 
aircraft  and  of  the  radar  parameters.  The  dissipation  factor,  as  used  in  tur¬ 
bulence  theory,  represents  the  kinetic  energy  converted  to  heat  per  unit  mass 
per  unit  time.  This  conversion  to  heat  occurs  at  the  end  of  the  sequence 
where  the  large  eddies  progressively  decay  into  smaller  and  smaller  eddies. 
The  conversion  from  kinetic  to  thermal  energy  occurs  on  the  scale  of  milli¬ 
meters.  MacCready^  has  related  e  to  the  response  of  a  particular  aircraft 
as  well  as  to  a  universal  turbulence  intensity  scale. 

♦This  work  is  an  extension  of  the  work  reported  in  Reference  1  and  is  noted  in 
Reference  2. 


The  quantity  e  can  be  determined  directly  by  the  aircraft  by  simply 

measuring  the  fluctuations  In  Instantaneous  airspeed.  Kolmogorov*  has  shown 
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that  the  structure  function  of  the  Instantaneous  airspeed  Dy  “  C(er)  ,  where 
Dv  is  defined  as  the  average  of  the  square  of  the  difference  of  two  successive 
airspeeds  v  measured  a  distance  r  apart.  It  Is  assumed  that  homogeneous 
Isotropic  turbulence  exists. 


Dy  -  E|(v2  -Vj)2}  -  C(er)2/3 


The  universal  constant  C  has  been  experimentally  determined  to  be  1.77  +  .08. 
It  can  also  be  shown  (see  Appendix  A)  that  the  acceleration*  an  aircraft 
experiences  Is  proportional  to  i.e., 


a r- 


KC1/2  pc1/3  1 
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Here  p  Is  the  air  density,  v  is  the  air-craft  speed,  m  Is  the  aircraft  mass 
and  K  is  a  constant  of  proportionality  unique  for  a  particular  aircraft.  It 
is  obvious  from  the  above  relations  that  the  quantity  rather  than  e 

Itself  is  the  more  useful  quantity. 
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The  ground-based  radar  can  also  determine  e  from  measurements  of  the 
width  of  the  rainstorm  radar  spectrum.  If  It  Is  assumed  as  before  that  the 
turbulence  is  isotropic  and  homogeneous,  then  the  following  relation  between 
et/3  and  the  spectral  width  ou  holds: 5 


r(l)o(ca)2/3F(^,  I;  1; 

3  3  2  2 


where  a  >  b 


(11 


♦Actually  the  square  root  of  the  acceleration  structure  function 


The  parameters  above  are  defined  as  follows: 

r(— )  "  gamma  function  (1.3541 1-—) 

‘  3 

a  -  universal  constant  (1.35  +  .06) 
a  -  radar  half  beamwldth  (cm) 
b  -  radar  half  pulse  length  (cm) 

F  -  Gaussian  hypergeometrld 

function 

2  —3 

e  *  dissipation  factor  (cm  sec  ) 

It  is  shown  in  Appendix  B  that  equation  [1]  can  be  closely  approximated  by 

o,.2  -  r(— )o( ea)2/3  -  1.828(ea)2/3  [21 

u  3 

Notice  that  Eqs.  [1]  and  [2]  do  not  contain  an  unknown  constant  of  pro¬ 
portionality.  It  is  for  this  reason  that  the  initial  emphasis  has  been  to 
measure  turbulence  in  terms  of  the  spectral  width  rather  than  some  other  pa¬ 
rameter  such  as  velocity  gradient  where  an  equivalent  theoretical  relationship 
is  not  available. 

3.0  EQUIPMENT  DESCRIPTION 

The  basic  experiment  consists  of  flying  an  instrumented  aircraft  through- 

a  turbulent  rainstorm  while,  at  the  same  time,  the  ground-based  lnstrumenta- 
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tlon  radar  attempts  to  estimate  the  amount  of  turbulence  (c  )  present  via  an 
appropriate  algorithm.  One  of  the  objects  of  this  program  is  to  determine  the 
best  radar  operating  parameters  and  algorithms  to  use.  A  brief  description  of 
the  aircraft  and  radar  instrumentation  follows. 

3.1  Aircraft  Equipment 

An  instrumented,  twin  engine,  turboprop  Grumman  Gulf stream,  operated  by 

the  FAA  Technical  Center  mi  used  to  penetrate  storm  cells  under  surveillance 
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by  the  ground-based  radar.  Within  the  storm  cells  the  aircraft  sensed  e  by 
measuring  the  pitot  differential  pressure  referenced  to  static  pressure  (IAS); 
outside  total  temperature,  and  absolute  pressure  (barometric  altitude)  were 
also  measured.  In  addition,  the  aircraft  carried  a  vertical  accelerometer 
mounted  at  the  center  of  gravity  of  the  aircraft  to  provide  a  measure  of 
relative  turbulence  level,  an  Inertial  navigation  system  to  provide  the 
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aircraft  position,  and  a  time-of-day  dock.  All  of  the  sensor  outputs  are 
digitally  recorded  on  tape  with  a  1  second  update  rate  with  the  exception  of 
the  pitot  pressure  and  the  accele rone ter  outputs  which  are  saapled  100  and  20 
times  per  second  respectively. 

3.?  Ground  Equipment 

The  ground-based  facilities  included,  in  addition  to  the  instrumentation 
radar,  a  radar  beacon  Interrogator  and  a  station-keeping  radar  which  were 
used  respectively  for  the  control  of  the  aircraft  and  for  the  location  of 
storms.  The  instrumentation  radar  consisted  of  an  S-band  ASR-8  coherent 
transmitter  with  its  front  end  coupled  to  a  15-foot  parabolic  antenna  mounted 
on  a  digitally  controlled  pedestal.  A  custom-built  receiver  consisting  of  a 
linear  IF  strip  drives  a  pair  of  quadrature  detectors  which  in  turn  drive 
10-bit  A/D  converters.  The  output  of  the  A/D  converters  fills  a  fast  2S6K  by 
20-bit  buffer  during  a  sector  scan  of  the  antenna.  The  data  is  then  read  out 
on  tape  via  a  NOVA  computer.  The  antenna  has  a  1.6-degree  one-way  beamwidth 
and  is  capable  of  being  pointed  in  azimuth  or  elevation  either  manually  or  by 
computer  control. 

f 

An  ATCBI-4  beacon  radar  is  used  to  locate  the  aircraft.  Its  antenna  has^ 
a  2.4-degree  azimuthal  beamwidth,  is  mounted  on  the  station  keeping  radar 
antenna,  and  is  located  30  meters  from  the  instrumentation  radar. 

The  stat ion-keeping  radar  is  an  FAA  S-band  ASR-7  terminal  radar  connect¬ 
ed  to  a  ASR-5  fan-beam  antenna.  Its  output  is  presented  on  a  PPI  to  locate 
storms.  The  instrumentation  radar  cannot  be  easily  used  for  this  purpose. 

A  Data  Entry  2nd  Display  System  (DEDS)  presents  both  analog  video  and 
alphanumeric  characters  on  a  PPI-llke  display.  Normal  video  from  the  station 
keeping  radar  and  the  beacon  returns  are  presented  as  analog  video.  NOVA- 
processed  data  such  as  dBz  level,  mean  radial  velocity,  turbulence  levels 
(e  )  as  well  as  the  most  recent  minute's  worth  of  bsacon  returns  are  super¬ 
imposed  as  alphanumerlcs.  This  allows  the  operstor  to  determine  in  real-time 
whether  a  particular  area  is  hazardous  and  should  be  avoided,  or  whether  it 
should  be  penetrated  by  the  instrumented  aircraft  to  gather  data. 
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4.0  CALIBRATION  AND  ERROR  ASSESSMENT 


Sources  of  error  in  the  experiments  requiring  special  attention  were: 

a.  Beacon  Antenna  Pointing  (AZ) 

b.  Instrumentation  Antenna  Pointing  (AZ  +  EL) 

c.  Instrumentation  Radar  and  Antenna  Gain 

d.  Beacon  Range  Error 

e.  Synchronization  of  Aircraft  and  Radar  Clocks 

Beacon  antenna  pointing  was  calibrated  by  observing  a  transponder  on  the 
ground  whose  position  is  known  relative  to  the  radar.  The  Instrumentation 
Antenna  pointing  was  calibrated  by  directing  the  antenna  to  an  accurately 
calculated  position  of  the  sun.  By  noting  how  far  the  antenna  has  to  be  mov¬ 
ed  tp  be  precisely  on  boresight  with  the  sun  it  was  found  that  the  error  was 
small,  less  than  0.1  degree  over  a  period  of  approximately  six  hours.  The 
Instrumentation  Radar  antenna  gain  was  determined  by  measuring  the  solar  flux 
level  when  observing  the  sun  and  comparing  it  with  the  published  values  for 
that  time  and  day.  The  radar  receiver  gain  was  determined  by  injecting  a 
calibrated  noise  source  into  the  radar  waveguide  directional  coupler.  A  fix¬ 
ed  beacon  range  error  of  0.26  nml  was  found  when  comparing  simultaneous  ob¬ 
servation  of  the  aircraft  by  the  instrumented  radar  and  the  beacon.  The  er¬ 
ror  is  attributed  to  the  differences  in  delay  between  the  ground  and  aircraft 
beacon  transponders  and  has  been  removed.  It  is  obvious  that  the  ground  in¬ 
strumentation  computer  clock  must  be  carefully  synchronized  by  the  aircraft 
Instrumentation  clock.  This  is  done  before  takeoff  and  rechecked  several 
times  during  the  flight  by  voice  communication.  The  timing  should  be  better 
than  one  second  accuracy  and  should  not  cause  any  appreciable  error  lr  uie 
overall  results. 

5.0  DATA  REDUCTION  PROCEDURES 

5.1  Beacon  Algorithm 

The  purpose  of  the  beacop  data  reduction  algorithm  is  to  accurately  lo¬ 
cate  the  instrumented  aircraft  relative  to  the  storm.  However,  when  the  bea¬ 
con  data  was  taken,  it  was  found  that  some  of  the  data  was  contaminated  by 
spurious  returns  from  other  aircraft  having  the  same  code,  from  the  ground 


test  transponder,  and  from  undetected  errors  in  the  NOVA  program.  The  errors 
have  been  corrected,  however,  since  the  data  is  valuable.  A  successful  ef¬ 
fort  has  been  made  to  retrieve  the  contaminated  data. 

Fig.  la  shows  a  plot  of  the  raw  beacon  returns  from  the  July  3  flight. 
Returns  from  the  ground  have  already  been  removed.  The  wide  scatter  is  pri¬ 
marily  due  to  the  old  error  in  the  NOVA  program  that  considers  returns  in  ad¬ 
jacent  range  gates  as  separate  targets.  As  a  result  of  this,  a  corrective 
program  was  written  to  cluster  the  returns  on  the  basis  that  any  returns 
within  a  certain  area  and  time  box  are  the  same.  The  results  are  shown  in 
Fig.  lb.  Any  missing  returns  are  then  filled  in  using  a  fourth  order  least 
square  fit.  Missing  points  have  been  caused  by  propagation  effects  such  as 
the  orientation  of  the  aircraft  or  flying  in  the  cone  of  silence  of  the  bea¬ 
con  radar.  Other  causes  of  missing  returns  have  been  an  error  in  the  NOVA 
program  (which  has  been  since  corrected)  and  possible  mlxups  in  the  settings 
of  the  data  collection  equipment.  Filling  in  the  missing  points  is  needed  to 
allow  proper  smoothing  in  those  areas  adjacent  to  the  missing  points. 
Fig.  lc  shows  the  missing  points  replaced  while  Fig.  Id  has  been  smoothed  and 
interpolated.  The  smoothing  and  Interpolation  is  accomplished  using  a  Fast 
Fourier  Transform  (FFT)  technique.  The  x  and  y  components  of  the  beacon 
track  are  Fourier-transformed,  zeros  are  added  to  the  high  frequency  end  so 
that  the  total  number  of  FFT  points  are  increased  by  four,  the  transform  is 
then  multiplied  by  a  Gaussian  weighting  function,  the  inverse  transform  is 
then  taken,  resulting  in  a  smoothed  beacon  track  with  four  times  as  many 
points.  Rather  than  having  a  beacon  file  with  an  update  of  4.8  seconds,  it 
is  now  every  1.2  seconds  and  will  allow  a  smoother  merging  with  the  instru¬ 
mentation  radar  data. 

5.2  Instrumentation  Radar  Algorithms 

The  instrumentation  radar  records  on  tape  the  unprocessed  radar  data  as 
an  in-phase  and  quadrature  component  for  each  of  the  range-azimuth  cells  in  a 
selected  sector.  The  tapes  are  then  sent  to  Lincoln  Laboratory  for  proces¬ 
sing  in  the  following  manner.  The  raw  data,  which  represent  the  returns  from 
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Fig.  Id.  Smoothed  and  Interpolated. 
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Che  rain  droplets,  were  recorded  range  sequentially.  However,  the  processing 
requires  that  the  data  be  first  reordered  in  azimuth.  A  maximum  entropy 
estimation  is  then  wade  of  the  autocorrelation  lags.  The  lags  are  used  in  the 
computation  of  moments.  The  moments,  in  turn,  provide,  by  using  specific 
pulse  pair  algorithms,  the  intensity  (dBz),  the  mean  radial  velocity,  and  the 
spectral  width  (in  e  units).  The  maximum  entropy  estimation  of  the  lags  is 
used  because  it  provides  a  better  estimate  than  a  straightforward  calcula¬ 
tion.  This  is  because  the  entropy  calculation  requires  that  the  estimated 
correlation  lags  be  consistent  with  one  another  (this  Insures  that  the  spectra 
the  lags  represent  are  positive  for  all  frequencies).  A  set  of  autocorrela¬ 
tion  lags  are  computed  for  each  radar  cell.  Typically  each  cell  is  made  up  of 
204  returns  (I  and  Q  components).  The  expressions  for  converting  the  correla¬ 
tion  lags  to  moments  are  as  follows: 

Intensity  or  0-th  moment: 


dBz  -  10  log(k|RQ|) 
Mean  radial  velocity  or  1st  moment: 


X  -1  x.<Rl> 
v  • -  tan  l(— — — ) 

4*T  Re(Rj) 


Velocity  variance  or  2nd  moment: 


Kl 

24*  V  Kl 


In 


[3a] 


[3b] 


[3c] 


where  Rq  Is  the  n'th  autocorrelation  lag,  \  the  radar  wavelength,  T  is  the 
Interpulse  period,  and  k  a  constant  determined  from  the  radar  parameters.  The 
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above  algorithms  are  subject  to  a  minimum  slgnal-to-noise  value.  In  particu¬ 
lar,  the  2nd  moment  estimator  requires  at  least  a  3  dB  S/N.  The  S/N  estimator 
used  is 


R  R1/3 
S  o  2 

_  -  l/< -  -  1). 

N  _  4/3 

K1 

The  2nd  moment  and  S/N  estimators  were  derived  assuming  the  spectra  is 
Gausslan-shaped.  Monte  Carlo  tests  of  these  estimators  show  that  the  algo¬ 
rithms  are  valid  for  a  time  series  that  has  a  Gausslan-shaped  spectrum. 
Homogeneous  isotropic  turbulence  theoretically  is  expected  to  produce  a 
Gausslan-shaped  radar  spectrum. 

The  computations  produce  arrays  of  numbers  that  represent  the  moments  at 
relatively  coarse  discrete  uniform  positions  over  the  recorded  sector  area. 
In  the  azimuth  dimension  they  are  a  beamwldth  apart  while  in  range  they  are 
spaced  112.5,  225  or  450  meters.  Consequently,  quantization  occurs  when  one 
attempts  to  assign  a  value  to  the  moments  along  the  aircraft  track.  Smoothing 
and  interpolating  these  discrete  points  is  accomplished  using  a  two-dimension¬ 
al  FFT  procedure.  Fig.  2  and  Fig.  3  are  contour  plots  that  were  processed  in 
this  manner.  Superimposed  are  the  aircraft  tracks.  In  general  the  tracks  ex¬ 
tend  in  time  halfway  to  the  previous  scan  time  to  halfway  to  the  next  scan 

1/3 

time.  Comparing  dBz  and  e  contours  for  the  same  scan,  it  is  clear  that  the 

e1^3  contours  are  much  more  structured.  The  intensity  plots  are  smooth  and 

1/3  1/3 

featureless  compared  to  the  widths  (e  ).  The  e  contours  on  the  right 
side  have  low  spectral  widths,  which  imply  low  turbulence  levels  while  the 
left  side  Indicates  that  there  are  areas  of  strong  turbulence.  In  contrast, 
the  intensity  contours  give  no  indication  of  heavy  turbulence  on  the  left 
side.  Comparisons  with  the  aircraft  data  will  be  discussed  later. 

5.3  Aircraft  Data  Reduction 

The  effect  of  turbulence  on  the  pitot-equipped  aircraft  is  sensed  by 
measuring  the  rms  acceleration  and  the  dissipation  factor.  The  acceleration 
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is  measured  directly.  The  dissipation  factor  e  is  found  from  the  pitot  dif¬ 
ferential  pressure  dp,  the  total  (stagnation)  temperature  T  and  the  absolute 
pressure  P^  (free  stream  pressure  or  barometric  altitude),  and  is  as  follows: 


1/3 


1/3  2/3 

<_£>  <ri>  ■T, 


RTt 

T 


Y+2  Y-l  2/3 

c  1/ 2 ( l  +  iE)  3y  C  Cl  +  iE)  Y  -  l]  p 


where  the  dp  structure  function  is  defined  as 

D  -  E{ (dp(t  +  t)  -  dp(-r))2}, 
flp 


[4] 


R  -  2.87  x  106(5!L)2K'1  gas  constant, 
sec 


Y  “  1.4 


ratio  of  specific 
heats. 


t  -  0.2-1. 0 
(sec) 


time  between 
successive  dp 
measurements , 


and 


C  »  1.77 


Kolmogorov's 

constant. 


Eq.  [4]  takes  into  account  the  aircraft  altitude,  airspeed  and  Mach  number. 

Both  the  structure  function  and  the  aircraft  rms  vertical  acceleration  are 

continuously  averaged  over  a  7.5  second  period  using  a  cosine  squared  weight- 

1/3 

lng.  Examples  of  time  plots  of  r  and  rms  acceleration  are  shown  in  the  up¬ 
per  and  middle  curves  of  Figs.  4  through  7. 
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6.0  COMPARISON  OF  AIRCRAFT  AND  RADAR  DATA 

1/3 

The  lower  curves  of  Fig.  4  and  5  are  the  values  of  e  measured  by  the 
radar  using  Eqs.  [2]  and  [3].  This  is  of  the  July  17,  1980  flight  and  will  be 
examined  at  length  here  because  it  includes  cases  where  the  aircraft  passed 
through  low-  and  high-turbulence  areas.  The  aircraft  made  many  passes  at  a 
5,000  ft.  altitude  through  a  moderately  violent  thunderstorm  in  the  Atlantic 
City  area.  Many  other  flights  have  been  examined  but  not  in  as  much  detail. 
Figs.  6  and  7  are  of  the  corresponding  reflectivity  in  dBz.  The  numbered  dots 
above  the  radar  data  (lower  curve)  represent  the  midpoint  time  when  a  particu¬ 
lar  radar  scan  (or  snapshot)  was  taken.  If  the  aircraft  is  in  a  radar  cell, 
it  will  appear  at  this  time.  The  smallest  value  of  reflectivity  that  can  be 

measured  with  the  radar  is  8  dBz.  However,  it  is  felt  that  in  order  to 

1/3 

trust  the  width  algorithms  (radar  e  )  the  reflectivity  should  be  greater 
than  14  dBz. 

The  time  interval  17:08  to  17:15  Fig.  4  exhibits  a  reasonable  correla- 

1/3 

tion  between  radar  and  aircraft  e  while  the  corresponding  reflectivity 

Fig.  6  does  not.  The  radar  peak  at  17:09  is  contaminated  by  the  aircraft. 

A  spectral  analysis  of  the  radar  return  at  this  time,  range  and  azimuth  show 

the  aircraft  to  be  present.  The  peak  at  scan  12  (17:10:30)  is  real  (no  air- 

1/3 

craft)  and  correlates  with  the  rms  acceleration  and  the  aircraft  e  .  A 

scatter  plot  has  been  generated  over  this  time  interval  and  is  shown  in 

Fig.  8a.  Here  the  aircraft  is  compared  to  the  radar  e1^.  A  reasonable 

correlation  is  evident.  The  corresponding  correlation  coefficient  is  0.81. 

The  correlation  with  reflectivity  (Fig.  8b)  is  0.57.  If  the  entire  flight  of 

1/3 

17  July  1980  is  used,  then  the  correlation  is  0.50  using  spectral  width  (e  ) 
and  0.29  using  reflectivity.  The  respective  scatter  plots  are  shown  in  Fig. 
8c  and  8d.  In  this  calculation  the  data  was  edited  to  remove  those  parts  that 
were  known  to  be  contaminated  by  the  aircraft.  For  example,  at  time  17:47, 
scan  39,  the  large  lobe  is  due  to  the  aircraft  in  a  range  cell.  See  Fig.  5. 
Fig.  9a  and  9b  show  the  interval  17:18  to  17:31.  The  width  correlation  is 
0.48  while  the  reflectivity  is  0.10.  This  data  includes  the  large  anomalous 
peak  at  17:30:10. 


( 


18 


t 


Fig.  4.  Comparison  of  radar  e1^3  (bottom  curve)  with  aircraft  e3^3 
(top  curve)  and  aircraft  rms  acceleration  (middle  curve).  (Flight  of 
17  July  1980). 
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Fig.  5.  Continuation  of  curves  of  Fig.  4. 
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Fig.  6.  Comparison  of  radar  reflectivity  (bottom  curve)  with 
(top  curve)  and  aircraft  rms  acceleration  (middle  curve).  (Flight  < 
July  1980). 


of  curves  of  Fig.  6. 


Fig.  7 .  Continuation 
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1/3  1/3 

Fig.  8.  (a-d)  Scatter  plots  of  radar  e  and  reflectivity  vs  aircraft  e 

(Flight  of  17  July  1980).  (a)  Correlation  coefficient  ■  0.81,  (b)  correlation 

coefficient  ■  0.57,  (c)  correlation  coefficient  »  0.50,  (d)  correlation  co¬ 
efficient  “  0.29. 
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1/3 

The  large  e  peak  at  17:03:14  has  a  double  maximum  spectrum  (Fig. 
10a).  The  pulse  pair  algorithm  Eq.[3c]  will  fail  under  these  conditions  and 
will  produce  too  large  a  value.  The  spectra,  in  general,  have  various 
shapes.  The  pulse  pair  algorithm  is  designed  specifically  for  Gaussian 
shapes,  however,  sawtooth,  flat-top,  multimaxima  as  well  as  Gaussian  shapes 
are  observed.  Consequently  it  is  not  surprising  that  perfect  correlations  are 
not  observed.  It  is  reasonable  that  an  algorithm  that  takes  into  account 
the  various  shapes  should  be  used.  Many  of  the  excessive  values  of  can 

be  ascribed  to  either  multi-modes  or  to  flat-top  spectra.  The  anomalous  peak 
at  17:30:13  has  a  wide  flat-top  spectrum  (Fig.  10b). 

7.0  CONCLUSIONS,  AND  FUTURE  PLANS 

1/3 

The  spectral  width  (e  '  )  is  superior  to  the  reflectivity  in  determining 

the  presence  of  turbulence.  The  reflectivity  is  featureless  compared  to 

1/3 

spectral  width.  The  correlation  between  aircraft  and  radar  e  is  not  per¬ 
l/3 

feet  and  in  general  the  radar  e  is  higher  and  on  occasion  may  be  extreme. 
1/3 

The  extreme  radar  e  values  appear  to  be  caused  by  flat-top  or  dual-maxima 
spectra. 

Work  is  underway  to  design  new  algorithms  that  will  take  into  account  the 
spectral  shape.  Since  the  raw  radar  data  has  been  recorded  on  tape,  various 
algorithms  are  being  tested  to  determine  the  one  which  results  in  the  best 
correlation  with  the  aircraft  data.  Comparisons  between  flights  have  to  be 
made  and  an  overall  correlation  coefficient  has  to  be  generated. 
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1/3 

Fig.  9.  (a-b)  Scatter  plots  of  radar  e  .  (Flight  of  17  July  1980). 
(a)  Correlation  coefficient  *  0.48,  (b)  correlation  coefficient  ■  C  •  >7. 
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APPENDIX  k 

AIRCRAFT  ACCELERATION  STRUCTURE  FUNCTION 


The  lift  of  an  aircraft  is  given  by 


1  2 
L  -  _  C  v  S 
2  L 


where 


C^  ■  coefficient  of  lift, 
p  *  air  density, 
v  ■  aircraft  true  airspeed, 

S  «  wing  area. 

The  coefficient  of  lift,  C^,  is  a  linear  function  of  the  angle  of 
attack,  over  the  range  of  interest.  Therefore 

1  u  2 

L  -  _  (_  C  )pv  S 

2  v  L 

a 

where  CL  is  the  lift  curve  slope  and  u  is  the  vertical  component  of 
a 

airspeed.  Thus  ^ 

L  ■  _  C  pSuv 
2  L 

a 


Taking  uifferences  we  find 


AL  “  _  C  pS [uAv  +  vAu] 
2  L 


Since  u«v  and  (Ap)^  a  (Av)^  we  find 


AL  “  —  C  pS[vAu) 
2  L 


and  since  the  force  AL  *  mAa  where  m  is  the  aircraft  mass,  we  obtain  for  the 
acceleration  difference 

1  pS 

Aa  ■  —  C  —  vAu 

2  L  m 

a 


The  structure  function  is  defined  as  the  average  of  the  square  of  the 
differences 

,  2,  1  pS  2 

D  -  E{ (Aa)  }  -  (-  C  _  v)  D 
a  1  2  L  m  u 

a 

Kolmogorov's  hypothesis  states  that  for  lateral  velocities 


Du  «  4/3C(er)2/3 

where  C  -  universal  constant  (1.77) 
e  -  dissipation  factor 

r  ■  separation  distance  when  measuring  the  two  velocities 
However,  r  ■  vt  where  t  is  the  time  difference  between  successive  acceleration 
measurements.  Therefore,  we  find  for  the  square  root  of  the  acceleration 
structure  function 


/F 


1  1/2  1/3  v 

—  CSC  p  e  — 


4/3 


✓? 


4/3 

1/2  1/3  v 

-  KC  pe  - 

m 

For  the  Gulfstream 

C^  -  5  (radians  *) 
a 

S  ■  610  square  feet 

m  *  933  slugs 


*S.  Panchev,  "Random  Functions  and  Turbulence,"  Pergaaon  Press  (1971)  p.  151. 
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SOME  BASIC  RELATIONS  CONCERNING  THE  RADAR  MEASUREMENT  OF  AIR  TURBULENCE 


INTRODUCTION 

The  use  of  a  radar  to  measure  the  level  of  air  turbulence  in  a  rain  storm 
is  a  relatively  new  concept.  Consequently,  the  purpose  of  this  note  is  to  pre¬ 
sent  the  connections  between  what  a  radar  measures  and  the  physical  properties  of 
the  turbulent  atmosphere.  In  the  first  section  we  will  show  that  the  radar 
spectrum  represents  the  velocity  distribution  of  the  rain  drops  and  appears  to 
have  a  Gaussian  shape.  In  the  second  section  we  will  show  how  the  radar  spectrum 
width  (variance)  is  related  to  the  turbulence  intensity  and  to  the  radar  cell 
size. 

SHAPE  OF  THE  RADAR  SPECTRUM 

Consider  a  volume  of  scatterers  such  as  rain  drops  that  are  being  carried 
by  a  homogeneous,  turbulent  atmosphere.  Let  us  look  first  at  one  of  the  rain 
drops.  The  radar  return,  in  voltage,  at  time  "t"  and  "t  +  t"  is  proportional 
to: 


where  x  is  a  random  variable  representing  the  position  of  the  drop  relative  to 
the  radar  and  X  is  the  radar  wavelength.  The  correlation  function  of  the  return 
is  by  definition 


where 


is  a  random  variable  representing  the  average  velocity  over  the  interval  t.  "t" 
drops  out  because  the  statistics  are  assumed  to  be  stationary.  Evaluating  (2) 
we  find 
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P(T)  "  /°°s<w)  elwT  da)  (4) 


where  S(u>)  is  the  spectrum  observed  by  the  redsr.  Equating  (3)  and  (4)  we  find 

,4ir 

.  i— r  UT 

S(uj)  e  WT  dw  •  2ir  p(u)  e  du. 


We  see  that 


4n  .  dm  4tr 

w  “  T  u  and  d^  "T 


and  find  that 


S(u)du  ■  2irp(u)du  and  S(u>)  «  j  p(u). 


(4a) 


(5) 


Thus,  we  see  that  the  spectrum  is  a  replica  of  the  probability  density  of  the 
velocity  averaged  over  the  interval  T.  Because  the  spectra  of  Independent  pro¬ 
cesses  can  be  added,  Eq.  (5)  is  valid  for  many  raindrops,  if  one  assumes  all 
the  droplet  cross  sections  are  the  same. 


It  is  not  difficult  to  generalize  Eq.  (5)  for  varying  droplet  sizes.  The 
return  at  time  t  for  N  droplets  Is 


N 

I 

n-1 


An  e 


*n(t> 


where  Aq  is  the  amplitude  of  the  nth  droplet.  The  covariance  function  Is 


N 

R(t)  -  E  {  Z  A.  e 
1-1  1 


141TX,  (t) 

F1 


N 

l 

J-l 


i4*x  . (t*T) 

.  X  *  , 

Aj  e  } 


Remembering  that  if  two  random  variables  (rv's)  are  independent,  the  expectation 
of  the  product  of  two  rv's  is  the  product  of  the  expectations,  we  find 

.4tt 

9  1_T  UT 

R(t)  -  N  E  {A  }  E  {e  A  } 


and 


R(0)  -  N  E  {A  } 


Since  x^  and  x^  are  independent,  the  cross  terms  drop  out  because 


i-^Y  (x,~x.|)  i^Y  *4  7  *4 

E{e  3  }  •  E  (e  }E{e  3 }  ■  0  when  j  j*  i . 


.4tt 
i-r  x. 


(6) 


Thus 


m 

/  S(a>)  e 


lUT 


du)  “  p<t>  *  f(oj- "  E  {e  X  } 


<*>  i^Y  ut 

/  p(u)  e  du 


or  as  before 

S(w)  -  y  p(u).  (5) 

Thus,  the  "radar"  spectrum  is  a  replica  of  the  velocity  probability  density. 
There  appears  to  be  some  theoretical  Justification^  for  the  velocity  distribution 
to  be  Gaussian.  Experimentally,  in  homogeneous  turbulence,  the  velocity  dis¬ 
tribution  is  symmetric  and  tends  to  be  nearly  Gaussian^1’2’33.  It  should  be  pointed 
out  that  in  a  related  process,  the  motion  of  a  molecule  in  a  gas,  that  the  dis¬ 
tribution  of  one  component  of  the  velocity  is  Gaussian  and  is  known  as  a  Maxwellian 
distribution. 

VARIANCE  OF  THE  RADAR  SPECTRUM 

Me  have  just  established  that  the  radar  spectrum  is  a  copy  of  the  velocity 

distribution.  It  follows  from  Eq.  (4a)  that  the  variance  of  the  radar  spectrum 

2  2 
o  is  related  to  the  velocity  variance  a  as  follows: 
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Since  E{u2}  is  not  a  function  of  r  and  ufTj)  and  u(?2)  are  correlated,  we 

have 

°i|2  8  E  tu2}  -  //  Ajtfj)  A(?2)  Bu  &v  r2)  d?j  d?2  (?) 

where  B  (itj,  r2)  is  the  correlation  of  the  two  longitudinal  velocities  at 
points  and  Because  the  turbulence  is  assumed  homogeneous 


hi  (ri*  r2)  8  hi  (ri  -  r2>- 

Since  is  a  correlation,  it  can  be  defined  in  terms  of  a  3-dimensional  Fourier 
transform 

i  t  *(r.  -  r_) 

Bn  (ri  '  r25  8  Ihi^  e  d^* 

The  energy  spectrum  of  the  3  dimensional  velocity  field  *(£)  is  often 

confused  with  the  radar  Doppler  spectrum  S(w).  4>. 0  (k)  follows  an  inverse 

4  5  ** 

11/3  power  law  *  and  has  a  pole  at  the  origin  while  S(u)  is  gaussian-shaped. 

#(£)  is  the  velocity  spectrum  while  S(«)  represents  the  velocity  distribution. 
id  has  the  dimensions  of  radians/second  while  t  is  radians/meter.  The 
function  *ll&  is  known  when  the  turbulence  is  isotropic  and  in  the 
inertial  subrange.  We  should  notice  that 

E{u2}  -  Bu  (0)  -  /*utf)  d  1c 

consequently,  Eq.  (7)  becomes 

°u2  *  /  <«*  -  //*<?,)  *<?2)  /♦„<*>  •**  '(?i ' <•*  d?, 

*  /  *tt  (ft  [ 1  -  /*(?[)  e  1  drj  e  2  dr2  |  it  (8) 

Evaluating  the  inner  integral,  by  changing  to  rectangular  coordinates,  we  have 
+  i£*r. 

/A(rj)  e  drx  ■ 


777375  .  2 

(2ir)  b  a 


/// 


-[^4 ,  _4] 

e  »-  2  a  2  b  J 


ky  y  ♦ k*  *) 
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e 


dx  dy  dz 


2 


(2ir)3/2b  a2 


oo  -  X— -  +  i  k  y  00 - +  i  k  z 

I  «  2.2  3  „yJ.  2.2  - 


-  — ,  +  1  k  x 
*  |  e  2b2  *  dx 


dz 


-i  (V  t2  +  k2  .2  +  k2.2]. 


Therefore,  Eq.  (8)  becomes 


*»2  ■  /  *u(I>  t1  - 


r.  2  . 2  ^  .  2  2  .  2  2n 

-  Ux  b  +  k  a  +  kz  a  ] 


]* 


For  isotropic  turbulence  in  the  inertial  subrange 


4 


rW 


<{>  -  a  -  -i-  > 


where 


and 


4ir  k  k 

E(k)  -  a  e2/3  k'5/3 

k2  -  k  2  +  k  2  +  k  2 

x  y  * 

k  -  \t\ 

a  -  1.35  +  .06  (dimensionless)* 

e  *  the  dissipation  factor. 


Substituting  and  using  the  following  relations 

2  2  2 
k  z  -  Vl  cos*  $ 
x 

dk  »  k2  sln<(>  dk  d<*>  d6, 

we  convert  Eq.  (8)  to  polar  coordinates 


°u2  -  ?  f  A"5'3  ♦ 

u  **1'  0  0  0 


ae 


4ir 
2/3  it 


-kZ  (a2  sin2  ♦  +  b 


cos 


’] 


dddk  d<#> 


/"  .lo3  ♦  /"  k‘5/3  [l  -  1,2  ,ln'T  +  “•‘♦>1  did*. 

0  0  L  J 
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Using  the  following  solution  for  the  inner  integrel 


v— 1  ..  -y  jP)  dx  •  -  — —  y  "  p  r  .  v 


J  xv-A  (1  -  e 
0 


r  (  ^  > 


we  obtain 


,  2  -  |  ae2/3  T  (5/3)  /  sin  3  *  (a2  sin2  *  +  b2  cos2  <fr)1/3  d$ 

u  8  Q 


and  using  symmetry  we  find 


,  2  -  ~  a e2^3  T  (5/3)  /  sinJ$  (a2  sin2  $  +  b2  cos2  4>)3^3  d$ 

u  4  0 


Substituting  t  -  cos  we  find 


a  2  -  |  oe2/3  r  (5/3)  a2/3  (1  -  t)  (1  -  [l  -  t)1/3  t~1/2  dt. 

U  O  -  ^ 


Since  we  are  dealing  with  real  numbers,  the  Integral  is  valid  when 

h2 

0  l2!  -1, 

»* 

The  Integral  representation  of  the  Gaussian  hypergeometrlc  function  is  given  by 

/  l  .  .  \  _  r  (t)  f  .b"*l  / 1  ^  \  C-b-1  r  <  . _ \  “S  J 


F(a,  b;  c;  z)  -  - -  J  t6'1  (l-t)0-1*'1  (l-tz)_a  dt 

r  (b)  r  (c-b)  0 

where  R  (c)  >  R  (b)  >  0. 
e  e 

Therefore, 

au2  -  r  (  j  )  a(e  a)2/3  F  (-  j,  j;  1  -  ^  )  (1C 

a 

2  2 

where  b  <a 

2 

The  hypergeometric  series  converges  when  h  1  and  is  given  by  the  monotonlc 
series' 

*A11  the  terms  after  the  first  are  negative.  It  can  be  shown  that  there  are  no 
maxima  or  minima  when  0  <  h  <  1. 


1  1  5  .  h  h2_  _5 _  3 

F(-  3.  2J  2;  “  1  "  15  ~  105  1701  h  *  ' 


[h-  1  - 


Since 


F(-  -J,  y;  0)  -  1, 


F,_  i  I.  I.  n 

r'  39  2*  2* 


12  r  ( |  >  r(  j  ) 


^  .917 


and  0  <  h  <  1 


we  find  that 


.917  <  F(-  i  I;  |;  h)  <  1. 


Consequently,  the  hypergeometric  factor  in  Eq.  10  has  little  effect  and 
that  the  approximation 


°u  -  r  (f5  °  <e  8)2/3 

is  valid . 


For  the  case  when  the  beamwldth  is  smaller  than  the  pulse  length 


0 


1, 


we  can  substitute 

2 

t  ■  sin  $ 

into  Eq.  9  and  find  in  the  same  manner  as  before 


2 

au2  -  r<f)  a  (e  b)2/3  F(-  2;  1  -  ^  ).  (11) 

b 

The  series  F  (-  2;  j;  g)  - 

,  4  8  2  32  3  r-  -  i  *2  n 

1  *  15  8  ’  105  8  “  1701  8  “  *  *  ’  C8  1  “  b2  3 

converges  as  before  when  0  <  g  <  1,  is  also  monotonic,  and  Is  bounded  by 
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55 


.491  <  F(-  i  2;  j;  g)  <  1 


Consequently,  the  approximation 


r  (f) 


a(E  b) 


2/3 


Is  not  always  justified. 

COMMENTS 

In  a  typical  radar  the  pulse  length  is  usually  smaller  than  the  beamvldth 

and  therefore  Eq.  10  will  apply.  The  spectral  width,  a  ,  will  thus  be  pro- 

1/3  u 

portlonal  to  (e  a)  .  Since  the  beamwidth  is  proportional  to  the  radar 

range,  and  that  the  buffeting  of  an  aircraft  can  be  shown  to  be  proportional  to 

1/3 

we  find  by  solving  Eq.  (10)  for  e  that 

_ 1 _ 

[rC2/3>  a  F(-±,  §;  f;  1  -  4  ) 

a  J 

O 

where  R  is  the  range  to  the  turbulence  and  C  is  a  constant.  Thus  the  buffeting  will 
be  directly  proportional  to  the  spectral  width  as  seen  at  the  radar  and  inversely 
proportional  to  the  cube  root  of  the  range. 

The  use  of  a  Gaussian  function  to  approximate  the  range  dimension  of  the 
cell  is  reasonable  if  one  considers  the  effect  of  convolving  the  Impulse 
response  of  the  radar  receiver  with  the  transmitted  rectangular  pulse  shape. 

The  effective  pulse  shape  will  be  rounded  and  look  like  a  Gaussian  curve. 

This  will  be  particularly  true  if  the  receiver  bandpass  is  closely  matched  to 
the  transmitted  pulse.  The  use  of  a  Gaussian  curve  to  approximate  the  main- 
lobe  of  the  antenna  has  long  been  used  and  is  considered  to  be  quite  accurate. 

The  variance  derivation  assumes  that  (1)  the  turbulence  is  uniform  through¬ 
out  the  cell,  (2)  we  are  in  the  inertial  subrange,  (3)  the  raindrops  move  with 
the  wind,  and  (4)  there  is  no  wind-shear.  These  assumptions  are  valid  be¬ 
cause  of  the  following  reasons.  (1)  By  using  small  radar  cells,  homogeneity 


,1/3. 


1/3 
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becomes  more  likely.  This  is  accomplished  by  operating  at  shorter  radar 
ranges  and/or  using  narrower  antenna  beams.  (2)  The  inertial  subrange 
signifies  that  the  effects  of  inertia  override  the  effects  of  viscosity  and 
of  gravity.  The  inertial  subrange  scale  extends  from  about  a  millimeter  to 
about  a  kilometer.  Consequently  radar  cells  less  than  a  kilometer  long  can 
be  filled  with  homogeneous  turbulence.  (3)  The  horizontal  velocity  com¬ 
ponents  of  the  rain  are  expected  to  follow  the  wind  but  not  the  vertical 
components.  At  most  antenna  elevation  angles  of  interest  (i.e.,  small  eleva¬ 
tion  angles),  the  horizontal  components  will  dominate.  (4)  First-order  ef¬ 
fects  of  wind- shear  broadening  can  be  taken  into  account  by  measuring  the 
mean  velocity  of  adjacent  cells  and  subtracting  out  the  effect. 

In  the  above  derivation  of  the  radar  spectrum  variance,  the  approach  is 
similar  to  that  of  Frisch  and  Clifford**  with  the  following  differences.  It 
was  not  necessary  to  explicitly  take  the  3-dimensional  Fourier  transform 
Fp(k)  in  £  space  of  the  antenna  pattern.  The  universal  constant  a  (A  in  the 
Frisch  paper)  is  not  .47  but  1.35.  The  solution  of  the  case  where  the 
pulse  length  is  longer  than  the  beamwidth  is  in  error. 
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APPENDIX  C 


1/3 

AIRCRAFT  MEASUREMENT  OF  ei,J 


The  following  parameters  are  available  of  the  Gulfstream  aircraft  to 
measure  turbulence 

Pf  -  Static  or  free  stream  pressure 
Ap  -  Pitot  tube  differential  pressure 
T  -  Stagnaatlon  or  total  temperature 
Pf  is  taken  from  the  static  port  of  the  pitot  tube.  Ap  -  P  -P^, 

where  Pt  Is  the  total  or  stagnation  pressure.  Isentroplc  and  ideal  gas 
flow  conditions  are  assumed  permitting  the  following  relations  to  be  used. 


P 

t 


Y 


(1  + 


Y-l  2  y-l 

-  M  ) 
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(1) 


T 


t 


-  1  + 


Y-l 

2 


2 

M 


(2) 


2 


YRTf 


Here  M  is  the  Mach  number,  y  is  the  ratio  of  specific  heats,  R  is  the  gas 
constant,  and  u  is  the  free  stream  airspeed. 

Combining  (1)  and  (3)  we  have 


Y-l 


P 


f 


2  tRTf 


(*> 


38 


Even  under  conditions  of  heavy  turbulence,  the  fluctuations  of  ?t,  Tt, 
M  and  u  will  be  saall  coapared  to  their  respective  aean  values.  Consequently, 
we  can  take  differentials  of  (4),  substitute  differences  and  we  find 


The  structure  function  D  is  defined  as  the  expectation  of  the  differences 
of  successive  measurements  squared  or 

Du(t)  ^  e{ [u(t+t)  -u(t)]2} 


Therefore,  we  can  rewrite  (5)  as 


u2  D 


(RTf  ) 


(6) 


Kolmogorov's  Hypothesis  relates  l>u  to  r  as  follows: 

Du  -  C(eut)2/3  (7) 

where  t  is  the  time  between  measurements  and  C  is  a  universal  constant. 
Combining  (6)  and  (7)  we  find 


1/3 

e 


-1/2  -1/3  -4/3 
C  t  u  RT 

f 


(8) 
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However,  since  Che  temperature  sensor  does  not  measure  the  free  stream 
temperature  but  the  total  temperature  T£,  we  must  combine  (1),  (2) 

and  (3)  with  (8)  and  finally  find 


1/3 

e 


RT  1/3  2/3 

c_i>  cri)  'Vp 

T  2y 


Y-l  2/3 


c1/2(i  +  iE)  3Y  1(1  +  iE)  Y  -I]  P 


(9) 


where  Ap  is  taken  to  be  the  mean  of  the  differential  pitot  tube  pressure  over 
the  averaging  period  of  the  structure  function. 
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